package NormalityTest;
import java.util.Arrays;

import org.apache.commons.math.stat.descriptive.moment.StandardDeviation;

/**
 * @author Ray 
 * based on http://www.kevinotto.com/RSS/templates/Anderson-Darling Normality Test Calculator.xls
 * this class will calculate the Anderson-Darling test statistic that a sample of data is normal
 * Keep in mind the test assumes normality, and is looking for sufficient evidence to reject normality.										
 *	That is, a large p-value (often p > alpha = 0.05) would indicate normality.	
 */
public class ADNormality {
	
	//Number of Data points
	private int np;
	//Sample Mean
	private double sampleMean;
	//Sample Sigma, deviation
	private double sampleSigma;
	//p-value calculations
	private double pValues[];
	//AD test statistic
	private double ad;
	//AD* test statistic
	private double adstar;
	
	private double f1i[];
	private double minsF1i[];
	private double f2i[];
	private double si[];
	private double sumofSi=0;
	//final p-value to determine if it is greater than 0.0005
	private double finalpValues;
	private double criticalVal;
	public  double calculateCV(){
		criticalVal=0.752/(1+0.75/np+2.25/(np*np));
		return criticalVal;
	}
	public ADNormality(double[] sortedData){
		pValues=new double[4];
		f1i=new double[sortedData.length];
		minsF1i=new double[sortedData.length];
		si=new double[sortedData.length];
		sumofSi=0;
		init(sortedData);
	}
	private void init(double[] sortedData){
		this.np=sortedData.length;
		this.sampleMean=calSampleMean(sortedData);
		this.sampleSigma=calSampleSigma(sortedData);
		calF1iandMinsF1i(sortedData);
		f2i=minsF1i;
		calsi();		
		ad=-sumofSi/np-np;
		adstar=ad*(1+0.75/np+2.25/(np*np));
		calPvaules();
	}
	/**
	 * @return p-Values for ADtest
	 * if p>0.0005 it is normal distribution
	 */
	public double doAD(){		
		return getFinalpValues();
	}
	private void calPvaules(){
		for(int i=0;i<pValues.length;++i){
			pValues[i]=0;
		}
		if(adstar<13&&adstar>=0.6){
			pValues[0]=Math.exp(1.2937-5.709*adstar+0.0186*(adstar*adstar));
		}
		if(adstar<0.6&&adstar>=0.34){
			pValues[1]=Math.exp(0.9177-4.279*adstar-1.38*(adstar*adstar));
		}
		if(adstar<0.34&&adstar>=0.2){
			pValues[2]=1-Math.exp(-8.318+42.796*adstar-59.938*(adstar*adstar));
		}
		if(adstar<0.2){
			pValues[3]=1-Math.exp(-13.436+101.14*adstar-223.73*(adstar*adstar));
		}
		finalpValues=getmax(pValues);	
	}
	private double getmax(double[] pValues){
		double temp=0;
		for(double i:pValues){
			if(i>temp)
				temp=i;
		}
		return temp;
	}
	private void calsi(){
		for(int i=0;i<f1i.length;++i){
			si[i]=(double)(2*(i+1)-1)*((Math.log(f1i[i])+Math.log(f2i[i])));
			sumofSi+=si[i];		
		}
		
	}
	private void calF1iandMinsF1i(double[] sortedData){
		for(int i=0;i<sortedData.length;++i){
			double temp=(double)(sortedData[i]-sampleMean)/sampleSigma;
			f1i[i]=normDistv(temp, 0, 1);
			minsF1i[i]=1-f1i[i];
			
		}
		Arrays.sort(minsF1i);//sorting
	}
	/**
	 * @param sortedData
	 * @return average of array
	 */
	private double calSampleMean(double[] sortedData){

		double sum=0;
		for(int i=0; i<sortedData.length;++i){
			sum+=sortedData[i];
		}
		return sum/sortedData.length;
	}
	/**
	 * @param sortedData
	 * @return deviation of given a sorted array
	 */
	private double calSampleSigma(double[] sortedData){		
		StandardDeviation sd=new StandardDeviation();
		double ad=sd.evaluate(sortedData);
		return ad;		
	}
	/**
	 * @param sampleSigma
	 */
	public void setSampleSigma(double sampleSigma) {
		this.sampleSigma = sampleSigma;
	}
	/**
	 * @param X
	 * @param mean=0
	 * @param sigma=1
	 * @return
	 */
	private  double normDistv(double X, double mean, double sigma) {
		
		double res = 0;

		final double x = (X - mean) / sigma;

		if (x == 0) {
			res = 0.5;
		} else {
			final double oor2pi = 1 / (Math.sqrt(2 * 3.14159265358979323846));
			double t = 1 / (1 + 0.2316419 * Math.abs(x));
			t *= oor2pi
					* Math.exp(-0.5 * x * x)
					* (0.31938153 + t
							* (-0.356563782 + t
									* (1.781477937 + t
											* (-1.821255978 + t * 1.330274429))));

			if (x >= 0) {
				res = 1 - t;
			} else {
				res = t;
			}
		}
		return res;
	}
	
	public double getSampleMean() {
		return sampleMean;
	}
	public double getSampleSigma() {
		return sampleSigma;
	}
	public double getAd() {
		return ad;
	}
	public double getAdstar() {
		return adstar;
	}
	public double getFinalpValues() {
		return finalpValues;
	}
	public int getNp() {
		return np;
	}
	public double getSumofSi() {
		return sumofSi;
	}
	
}
